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In this paper, we analyze the temporal evolution of the age-dependent force of infection and 
incidence of rubella, after the introduction of a very specific vaccination programme in a previ- 
ously nonvaccinated population where rubella was in endemic steady state. We deduce an integral 
equation for the age-dependent force of infection, which depends on a number of parameters that 
can be estimated from the force of infection in steady state prior to the vaccination program. We 
present the results of our simulations, which are compared with observed data. We also examine 
the influence of contact patterns among members of a community on the age-dependent intensity 
of transmission of rubella and on the results of vaccination strategies. As an example of the the- 
ory proposed, we calculate the effects of vaccination strategies for four communities from Caieiras 
(Brazil), Huixquilucan (Mexico), Finland and the United Kingdom. The results for each community 
differ considerably according to the distinct intensity and pattern of transmission in the absence of 
vaccination. We conclude that this simple vaccination program is not very efficient (very slow) in 
the goal of eradicating the disease. This gives support to a mixed strategy, proposed by Massad et 
al, accepted and implemented by the government of the State of Sao Paulo, Brazil. 

PACS numbers: 87.10.+e, 87.19.Xx 



I. INTRODUCTION 

The control of directly transmitted, viral childhood in- 
fections, around the globe has been strongly dependent 
on vaccination, the most effective control tool developed 
so far [1]. There are several infections for which vac- 
cines exist. These are therefore candidates for eradica- 
tion. Some examples include polio, measles and rubella, 
just to mention a few. Vaccination strategies, however, 
have been each time more dependent on inferences based 
on quantitative models, which can, through simulation 
tools, yield distinct scenarios and possibilities. These 
simulation techniques, in turn, have been proved to be 
invaluable tools for helping health authorities to decide 
between competitive strategies of eradication or control 
of those infections. 

In previous publications 0, H[ , we applied mathemati- 
cal models to design and to evaluate the impact of vacci- 
nation against rubella in the state of Sao Paulo, Brazil. 
Rubella is a viral infection that causes a mild disease, 
but it is considered to be a public health problem due 
to the risk of fetal infection and subsequent congenital 
defects @, 0, H|. Therefore, the goal of rubella vaccina- 
tion is to prevent from the congenital rubella syndrome 
(CRS).Plotkin[| argues that, due to the high prevalence 
of rubella in some countries, only high vaccine coverage 
would avoid the increase of CRS. 

In this paper, we analyse the effects of different con- 
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tact patterns on vaccination strategies against rubella in 
some communities. We investigate a plausible form for 
the contact rate function including some constrains it 
must satisfy. We concentrate on the integral equation 
for the age-dependent force of infection — denned as the 
age-dependent number of new infections per capita, per 
unit time — , which relates the pattern of contacts among 
the members of a population with the prevalence of the 
disease, following a methodology developed elsewhere [f|. 
The basic idea is to examine the force of infection in 
steady state that results from a given vaccination strat- 
egy- 

We also turn our attention to the dynamics of the 
process, following the time development of the age- 
dependent force of infection when a vaccination strategy 
is started at a certain time t in a previously nonimmu- 
nized population. Some aspects of the age and time de- 
pendences in epidemic models have already been studied 
by some authors (e.g., @, [I])- 

This paper is organized as follows. In Sec. [Til we 
present the formalism used. We describe in detail how 
the contact rate function is related to the force of infec- 
tion, discuss some constrains it must satisfy, and propose 
a form for it. In Sec. Hill we describe the fitting proce- 
dures adopted to determine the values of the parameters 
of the contact rate function for different communities. 
In Sec. IIV Al we analyse the impact of specific vaccina- 
tion strategies against rubella using data from commu- 
nities from Caieiras, a Brazilian small town located in 
the neighbourhoods of Sao Paulo city (Azevedo Neto et 
al. Q), Huixquilucan, Mexico (Golubjatnikov et al. 0), 
Finland (Edmunds et al. [l(|) and the United Kingdom 
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(Farrington et al. [Ti|). It must be noted that the results 
from Brazil and Mexico are from nonvaccinated commu- 
nities, while the results from Finland and the United 
Kingdom are from nonvaccinated males in communitites 
that have partially vaccinated female population [ll], [l^] • 
The results from Sao Paulo will be compared with those 
previously reported by Massad et al. LL3| • This paper 
differs from those quoted above, in that the form of the 
contact rate function, representing the age related pat- 
tern of contacts, was studied more carefully. Also the 
relation between the vaccination rate v and the result- 
ing proportion of vaccinated people [see Eqs. pfi|) - (|T7)) ] 
was modified. In spite of this, as we shall see, the rec- 
ommended vaccination strategy was maintained, but the 
calculated effects of the vaccination strategies seem now 
to be more realistic. In Sec. IIVB1 we present simulations 
of the temporal evolution of the force of infection and, in 
Sec. IIV CI we compare our results to experimental re- 
sults. Finally, in Sec. [V] we summarize our results. 



II. MATHEMATICAL DEVELOPMENTS 
A. Temporal Evolution 

Let us assume a SIR model (Susceptible-Infected- 
Recovered). Let S(a,t)da, I(a,t)da, and R(a,t)da be, 
respectively, the number of susceptible, infected and non- 
susceptibles (including recovered and vaccinated) indi- 
viduals with ages between a and a + da at time t. We 
can write 



dS(a,t) dS{a,t) 



da 



dt 



dl(a,t) dl(a,t) 
da dt 
dR{a,t) dR(a,t) 



da 



dt 



~[X(a,t) + v(a,t) + fi]S(a,t) 
A(o,t)5(o,t)-(/i + 7)/(o,t) (1) 
u(a, t)S(a, t) + 7/(0, t) — /ii?(a, t) , 



where v(a, t) is the age and time-dependent rate of vacci- 
nation, 7 is the recovery rate, and fi is the mortality rate, 
assumed constant. This type of mortality rate (constant) 
is known as type-II mortality function. Another type of 
survival curve (type I) considers that all individuals sur- 
vive to exactly a certain age, and then die. Anderson 
and May [13] mention that, for both developed and de- 
veloping countries, the observed mortality function is in- 
termediate between type I and type II, although closer 
to type I for developed regions. 

The definition of the force of infection, as a function of 
age and time, is 



A(a,t) 



da p{a,a) 

N(a',t) 



(2) 



and N(a, t) = S(a, t)+I(a, t)+R(a, t) is the total number 
of individuals whose ages are between a and a+da at time 
t. In this equation, /3(a, a') is the so-called contact rate 
function. It is defined so that /3(a, a')dada' is the number 



of contacts a person with age between a and a+da makes 
with all persons with age between a' and a' + da' per unit 
time. Therefore, f](a,a') describes the contact patterns 
among the members of a population. 

Taking into account the three equations of system |T]), 
we can write, for iV(a, t), 



d d 



da dt 



+ —)N(a,t) = - l iN{a,t). 



(3) 



For simplicity, we consider that, at time t, the total 
population has size N. In other words, we have taken 
N(a,t) = N(a) = A(0)e _ '"\ for a given t. In this equi- 
librium situation, the mortality rate equals the natality 
rate, and we have N(Q) — fiN. 



1. Integral equation for \(a,t) 

Applying the method of the characteristics, as pro- 
posed by Trucco [l5[ (see also [l6j]) for solving the 
McKendrick-Von Foerster equation, we can solve the sys- 
tem of equations ([1]). 

Let s(a, t) and i(a, t) be the proportions of susceptible 
and infected individuals, among those with age a at time 
t, given by 



, . S(a,t) , . I(a,t) 
s( a , t) = ; ' , i(a,i)- 



N(a,t) 



N(a,t) 



(4) 



With these previous definitions, the first two equations 
of the partial differential equations system |T]) can also 
be written as follows: 



ds(a,t) ds(a,t) 

da dt 
di(a,t) di(a,t) 



da 



dt 



-[\(a,t)+v(a,t)]s(a,t) (5) 
A(a, i)s(a, t) — ji(a, t) . (6) 



The boundary conditions are such that, at age a — 0, 
for t > 0, we have s(0,t) — 1 and i(0,t) — 0. At time 
t = 0, for a > 0, we have that s(a, 0) and i(a, 0) are 
functions of age. In the calculations, the upper limit for 
the age is taken to be L = 60 yr. 

Considering the change of variables (as those proposed 
by Trucco [ijj) 



a — t 
t 



we have 



s(a,t) = s(£ + 77,77) = s'(£,?7) 



and similarly for i(a,t), \(a,t) and v(a,t). 

We also have that (d/da + d/dt) = d/drj. Thus, taking 
into account the above mentioned change of variables, 
Eq. §5§ reads 



d_ 

dt] 



In 77) = -[A' (£,77) + !/(£, 77)] 



(7) 
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whose generic solution can be written as 

hxa'(i,ri) = - r[\'(Z,x)+v'(Z,x)]dx + f(0 , (8) 
Jp 

where p and /(£), parameters related to the boundary 
conditions, are given by 



f(0 = Ina'foO), p = 
/(0 = Ina'^.-O, P=-£ 



(9) 
(10) 



for £ > and £ < 0, respectively. 

Then, for the cases in which £ > (a > t) or £ < 
(a < i), we have, respectively, the following solutions: 

£>0: s'($,ti) = s'($,0)e-fi [X 'G'* )+ ' /i t' x)]dx (11) 
£<0: = S '(£,-e)e"A [V( ^' ( ^ )lda: (l2) 

Rewriting the equations above in terms of a and t, we 
obtain 



[A(a — £ + x, x) 



(13) 



[A(a — t + x, x) 



s(a,t) = s(a — i,0)exp 

+ z/(a — t + a;, a;)] dx 

s(a,t) = s(0, t — a) exp 

+ z/(a — i + x, x)] dx 

= s(0, t — a) exp 

+ 1/(2, z — a + t)] dz 

for a > t and t > a, respectively. 
Equation © 



^ + ^ J i(o, t) + ji(a, t) = A(o, t)«(a, t) (15) 



can be rewritten, with the change of variables, as 



[X(z, x — a + t) 



(14) 



0,/ 



<9rj 

whose solution is 
■Ut ^ ~ f" 7ds 



i'{^ri)+li'{i,ri)=\ l {i,ri)s'{i,ri) , (16) 



(17) 



where q and g(£) depend on the boundary conditions: 

9(0 = i%0), q = (18) 



g(0 = i'(t,-0, q = -Z 



(19) 



for £ > and £ < 0, respectively. 

Equation (jTTJ) , in terms of a and t, is given by 

i(a,t) = [ dt' X(a-t + t',t')s(a-t,0)a(a,t,t') + 
Jo 

+ e~ 7t i(a - t,0), a>t, 



i(a, t) 



(20) 



da 1 A(a , a' — a + t)s(0, t — a)ip(a, t, a') + 



(21) 



+ e" 7Q i(0,i - a), a < t, 
where a(a,t,t') and tf)(a,t,af) are 

a(o, f , = e~ ^0' [A(a-t+T,r)+,(«-t+T,r)]dr e7 ( t '_i) (22) 

and 

V>(a < a') = e - -^" [*(*.«-<*+*)+»'(*.* ! -«+t)J<k e7 (a'-a) 

(23) 

The age and time-dependent force of infection [Eq. ^j] 
can also be written as 



X(a,t) = / da'/3(a,a')i(a',t) 
Jo 



(24) 



In the calculations, as already explained, the upper 
limit of the above integral is taken to be L — 60 yr. 
Thus, replacing solutions (|2H)l and (f2"Tj) for i(a,t) in the 
above definition, and considering that age is in the in- 
terval < a < L, the integral equation for the age and 
time-dependent force of infection is given by 



X(a,t) = 



i(t,L) 



da (3 (a, a) I da" X(a" , a" — a' + t) 



x i/)(a',t,a") +0(L-t) / da' p '(a, a') 



dt' X(a -t + t', t')s(a' - t, 0)a(a, t, t') 



+ e" 7 * i(a' -t,0) 



(25) 



where 8(L — t) is the Heaviside function. In the following 
section, we study the steady state of Eq. (1251) . 



2. Steady state behavior 

Let S(a)da be the number of susceptible individuals 
with age between a and a + da. The fraction of poten- 
tially infectious contacts they make with actually infec- 
tives aged between a' and a' + da' per unit time is 



S(a)daf3(a, a')da 



, Ha/) 
N(a') 



(26) 
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The total number of potentially infective contacts of 
susceptibles aged between a and a + da with infectives 
can be obtained by integrating Eq. (|26| in da' . Then, 
we obtain an expression for the age-dependent force of 
infection similar to Eq. ([24)) . 

Equation ([21) in the steady state condition gives 



i(o) = e^ a 



da'e 7Q 'A(a')s(0) 



x e 



(27) 



Substituting this expression in the definition of the age- 
dependent force of infection in steady state we have 



A(a) 



da'l3{ ai a') / da"e-^ a '- a "h{a") 



x e 



- jf dz[X(z)+y(z)} 



(28) 



The integral equation ([28)) always has A(a) = as solu- 
tion. According to Lopez and Coutinho [l7| . depending 
on the parameters of the integral equation, it may have 
another unique positive solution. 

Equation ([2"8"]) is the limit for large t of Eq. ([231) : 



A(a) = lim A(a, t) . 

t — >oo 



B. Contact Patterns 

1. Symmetry in the contact pattern 

As mentioned in the Introduction, one of our main diffi- 
culties is to choose a correct form for the contact function 
/3(<z, a'). In this section, we analyze a specific situation in 
which j3(a, a') has to satisfy a symmetry relation that 
restricts its form: if a person A has a contact with a 
person B, then B had a contact with A. In terms of 
transmission dynamics, it means that the total number 
of contacts a group C of infected individuals make with 
a group D of susceptibles equals the number of contacts 
group D had with group C. This symmetry is relevant 
when a direct, person-to-person contact is required for 
transmission. For instance, a direct contact is required 
for sexually transmitted diseases. It seems to be at least 
partially required for the transmission of directly trans- 
mitted childhood diseases such as rubella. 

The number of contacts the susceptibles with age be- 
tween a and a + da make with infectives with age between 
a' and a' + da' , in a time interval dt, is, as we have seen, 



S(a)da/3(a,a')da'j^dt 



(29) 



This number must be equal to the number of contacts 
the infectives with age between a' and a' + da' make with 



the susceptibles with age between a and a + da. This 
number is 



I(a')da'P(a', a)da^p-dt 



S(a) 

'N( a y 



(30) 



Thus, we must have 



S(a)P(a, a')^ = l(a')0(a', a)j& (31) 



P(a,a') _ (3{a',a) 
N(a') ~ N(a) 



(32) 



Since N(a) = A(0)e _ ^ a , we see that Eq. (J32J) is satis- 
fied if f3(a,a') has the form 



I3(a,a') = e^ a h{ ai a') 
where h(a, a') is symmetric, that is 
h(a, a') — h(a' , a) 



(33) 



(34) 



Equation ([33)) will be used in the following section to 
construct an analytical form for (3{a, a'). 



2. A form for the contact function f3(a, a') 

Let us consider that rubella is approximately trans- 
mitted by direct person-to-person contact. In this case, 
considering that children are stratified mainly by age in 
classrooms [3J|, it is reasonable to assume that contacts 
are more intense among children with the same age. It 
is then convenient to write h(a,a') as a product of two 
functions, 



h(a,a') = f(a,a')g(a,a') 



(35) 



The function f{a, a') represents the longitudinal profile of 
h(a, a') along the plane a = a! and g(a, a') the transversal 
profile related to the spread of h(a, a') to both sides of 
the plane a = a' . 

We have chosen the following positively skewed func- 
tion for f(a, a') 

f(a,a') = b 1 (a + a')e' b ^ a+a '^ , (36) 
and a Gaussian-like function for g(a,a') 

s -(a~a') 2 /* 2 



g(a,a 



(37) 



where a = a (a, a') is related to the width of the 
Gaussian-like distribution to the sides of a = a'. Consid- 
ering a linear spread 

a (a, a') = b 3 + b 4 (a + a') , (38) 

we obtain 

h(a, a') = b 1 {a + a')e- b ^ a+a '^ e - {a - a ' )2 V^+M^+a')] 2 ^ 

(39) 
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where bi, 62, 63 e 64 are the parameters to be determined. 

Thus, taking into account Eq. ([33]) . we have, for the 
contact function fl(a,a'), 

(3(a, a') = 61 (a + a ')e~ b2( - a+a ' ) e~ {a ~ a " >2 /[b3+bi{a+a " >]2 . 

(40) 

Other functions could be chosen for h(a,a'), as those 
proposed in Coutinho et al. 6] and Massad et al. 



C. The relationship between vaccination rate and 
vaccine coverage 

For our next simulations, we need to define what we 
mean by vaccination routine. In a nutshell, we take 



u(a, t) — i>6{a — ao)9(ai — a)9(t — to) 



(41) 



which has the following interpretation: after time to 
years, children are vaccinated at a constant rate of v 
children per unit of time when their ages are between 
ao and a\. In practice, the government usually informs 
through the media that mothers should take their chil- 
dren to health centers to receive the shots. The response 
of parents to the government advertisement results in a 
given v. Enthusiastic response results in a high v. 
In steady state, Eq. ((Til) becomes 



(42) 



v(a) = 1/6 (a — ao)0(ai — a) 



We shall now calculate the relationship between v and 
resulting proportion of vaccine coverage, p. Let V(a)da 
be the number of vaccinated individuals with age be- 
tween a and a + da. Let N v (a)da be the number of non- 
vaccinated persons with age between a and a + da. We 
have 



dV(a) 
da 
dN v (a) 
da 



= v{a)N v (a) - fiV(a) (43) 
= -v(a)N v {a) - nN v (a) . (44) 



Of course, we have V(a) + N v (a) — N(a). 

Solving Eq. (|43[) using the form of u(a) given by Eq. 
(|42l). we have 



!0 , a < a 

7V(0)e^ Q [l - e -Ko-<»o)] ; a < a < a x 
7V(0)e-' ia [l - e -"(ai-ao)] j a > ai 

The proportion p of vaccine coverage is defined as 

V ( a l) _ 1 „- I /(ai-o ) 



P 



1 - e 



N( ai ) 

The inverse relation between v and p is 
ln(l-p) 



(45) 
(46) 



a - a\ 



(47) 
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FIG. 1: Seroprevalence data and corresponding fitted curves 
for communities from Brazil, Mexico, Finland and the UK. 
The data for Caieiras (Brazil), Huixquilucan (Mexico), Fin- 
land and the UK were taken, respectively, from the works of 
Azevedo Neto et al. [Hi Golubjatnikov et al. Edmunds et 
al. [10] and Farrington et al. [11J. 



III. FITTING THE MODEL TO THE DATA 

Data consisted in seroprevalence studies carried out in 
communities from Mexico, Brazil, Finland and the UK. 

Let S + (a)da be the proportion of seropositive individ- 
uals to rubella — whose serological tests were positive, 
indicating that they have already been infected — with 
ages between a and a + da. An estimate of the function 
S + (a) resulted from fitting the serological data to (see 
Ref. (HI) 



S + (a 



= l-exp|^[(A ;2a +l) e - fc2Q -l]| 



(48) 



where fe« (i — 1,2) are fitting parameters, estimated by 
the maximum likelihood technique for all the communi- 
ties except that from Finland, which was estimated by 
the least squares fitting technique. Figure [1] shows the 
results of the fitting functions for the four communities 
considered, and the fitting parameters are shown in Table 

m 

In our model, the seropositive individuals correspond 
to those who are either infected or nonsusceptibles (re- 
covered and vaccinated), i.e., the proportion of seroposi- 
tives, S + (a), is equivalent to 1 — s(a). The force of infec- 
tion in the absence of vaccination, Arj(a), was estimated 
from the seroprevalence data by the so-called catalytic 
approach (e.g., Ref. [l9(), according to 



X (a) 



dS+{a) 
da 



(l-S+(a))" 



(49) 
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The term catalytic arises from an analogy with chem- 
istry. In the dynamics of infectious diseases, an infected 
individual would act as a catalyst, infecting susceptible 
individuals. Equation (|4"5)) corresponds to Eq. ([5]) in 
the steady state for the susceptible individuals, in the 
absence of vaccination. 

Equation (14T))) , expressed in terms of Eq. ([38J) , results 

in 



Ao(a) = fciaexp [— fea] 



(50) 



The values of the parameters of the contact function 
/3(a,a') [Eq. ([40]) ] were calculated so that the resulting 
force of infection A (a), in the absence of vaccination, ob- 
tained by solving Eq. ([28]) iteratively, agreed with Ao (a) 
given by Eq. (|50|) . The parameters 7 and /i were taken, 
respectively, to be 26.0 yr -1 , corresponding to an infec- 
tious period of 2 weeks, and 0.017yr _1 , the inverse of a 
life expectancy of 60 yr. Those parameters were taken 
to be the same for all communities, for simplicity. The 
resulting parameters of the contact function /3(a, a') for 
each community considered arc shown in Table [IJ 

For Finland and the UK, we carried out simulations 
considering the two types of mortality functions de- 
scribed in Sec. Ill Al As the results were very similar, 
we discuss only those concerning type-II mortality rate. 



TABLE I: Fitting parameters (fci and fe) of the seropreva- 
lence function, and the parameters of the contact function, 
for each community considered. 



Community fci(yr 2 ) 


Myr- 1 ) 




b2(yr- 1 ) 




bi 


Brazil 


0.0456 


0.108 


0.658 


0.0468 


3.49 


0.341 


Mexico 


0.214 


0.255 


3.54 


0.116 


1.04 


0.416 


Finland 


0.0290 


0.1068 


0.587 


0.0608 


2.77 


0.398 


UK 


0.0833 


0.1804 


1.60 


0.0928 


1.747 


0.391 



The forces of infection [as given by Eq. (|50|) ] for the 
same communities are shown in Fig. [5J As can be noted, 
the curves have strikingly different shapes, reflecting dis- 
tinct contact patterns. As we shall see later, this has 
profound impact on the calculated efficacy of different 
vaccination strategies. 

From the force of infection, we can define the average 
age at which susceptibles acquire infection 



/ °° a\(a)s(a)da 
J °° X(a)s(a)da 



(51) 



We have taken the highest ages observed in the seroepi- 
demiological studies as the upper integration limits of 
the integrals of Eq. (f5Tj) . The calculated values for the 
communities studied are given in the Table HT1 below. 

The contact functions /?(a, a') [Eq. (|4"0"]) ] of the com- 
munities considered are shown in Fig. [3] as examples of 
the general shape obtained. The analysis of these contact 
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FIG. 2: Force of infection for the communities studied derived 
from equation (|50[) . 



TABLE II: Average age at infection for the four communities 
studied. 



Community 


a (yr) 


Caieiras, Brazil 


8.45 


Huixquilucan, Mexico 


3.96 


Finland 


10.6 


UK 


6.64 



functions suggests two distinct patterns. In Mexico and 
in the United Kingdom, the age distribution of contacts is 
concentrated at lower ages. In contrast, the communities 
of Caieiras and Finland show a broader range of contacts, 
spread over all ages. In addition, it can be noted that 
the density of contacts estimated for the communities of 
Mexico and Caieiras are roughly twice as high as in the 
United Kingdom and Finland, respectively. This may re- 
flect distinct social contexts between the developed and 
developing countries as well as the fact that data from 
developed communities are only for males in communi- 
ties which have partially vaccinated female populations 



IV. SIMULATION RESULTS 
A. Effects of specific vaccination strategies 

We now calculate the results of specific vaccination 
strategies in the above mentioned communities, choosing 
ao = 1 year and a\ = 2 years, ao = 7 years and a\ = 8 
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FIG. 3: Calculated values of the contact functions /3(a, a') and respective contour plots for: (a) Caieiras, Brazil; (b) Mexico; 
(c) Finland; and (d) the United Kingdom. 



years, and a = 14 years and a\ = 15 years for several 
values of v. For a given vaccination coverage proportion 
p, we determine v through Eq. (1471) . 

The simulated results of the vaccination strategies were 
obtained by solving Eq. using the values of the pa- 
rameters of f3(a,a') obtained in Sec. IIIII for the vaccina- 



tion strategies described above. 

The results for the communities of Brazil and Finland 
are shown, respectively, in Figures [4] and [5] The results 
for Mexico and UK are not shown in graphs, but we have 
discussed them below. Figure[3]shows the results of vacci- 
nation strategies applied to the community of Sao Paulo. 
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FIG. 4: Effects of different vaccination programs calculated 
for Caieiras, Brazil. Children are vaccinated (a) between 1 
and 2 yr , (b) between 7 and 8 yr, and (c) between 14 and 
15 yr. The numbers above the dashed lines indicate the cor- 
responding vaccine coverage and the solid lines correspond to 
the catalytic model. In graph (c), the dashed lines correspond, 
respectively, to 0.0, 0.4, 0.8 and 0.97 vaccine coverages. 



FIG. 5: Effects of different vaccination programs calculated 
for Finland. Children are vaccinated (a) between 1 and 2 
yr, (b) between 7 and 8 yr and (c) between 14 and 15 yr. 
The numbers over the dashed lines indicate the corresponding 
vaccine coverage and the solid lines correspond to the catalytic 
model. In graph (c), the dashed lines correspond, respectively, 
to 0.0, 0.4, 0.8 and 0.97 vaccine coverages. 



Figures 4(a)-4(c) represent the different age intervals of 
vaccination. It can be noted that 75% coverage in the 
age interval from 1 to 2 yr almost eliminates the disease, 
but the peak of infection is shifted to around 17 yr, and 
therefore it is displaced to the right as compared with 
the case of no vaccination. A coverage between 79% and 
80% eliminates the disease. In Fig. 4(b), it can be noted 
that 85% coverage in the age 7-8 yr almost eliminates the 
disease. In addition, the peak of infection occurs around 
8 yr, and therefore it is displaced to the left as compared 
with the case of no vaccination. A 90% coverage in this 
age interval eliminates the disease. Finally, Figure 4(c) 
shows that vaccination in the interval from 14 to 15 yr is 
almost useless, since 97% coverage has very little impact 
in the force of infection, and it is impossible to eliminate 
the disease, even if a 100% coverage is used. 

Figure [5] shows the results of vaccination strategies ap- 
plied to the community in Finland. Figures 5(a)-5(c) 
represent the different age intervals of vaccination. It 
can be noted that 60% coverage in the age interval 1-2 
yr [Fig. 5(a)] almost eliminates the disease, and the age 
of the peak of infection is not affected at all. A cover- 



age of 64% eliminates the disease. In Fig. 5(b), it can be 
noted that 70% coverage in the age interval 7-8 yr almost 
eliminates the disease, and again does not shift the age 
of the peak in the force of infection. Finally, Fig. 5(c) 
shows that vaccination in the interval from 14 to 15 yr is 
almost useless, since 97% coverage has very little impact 
on the force of infection, and it is impossible to eliminate 
the disease even if a 100% coverage is used. 

For the Huixquilucan community in Mexico a 74% cov- 
erage in the age interval 1-2 yr eliminates the disease. 
Even a 97% coverage in the age interval 7-8 yr is not 
able to eliminate the disease, and indeed causes very lit- 
tle effect on its force of infection. 

For the community in the UK, a 66% coverage in the 
age interval 1-2 yr eliminates the disease. Even a 97% 
coverage in the age interval 7-8 yr is not able to eliminate 
the disease. However, the peak of the force of infection 
curve shifts leftwards to around 5 yr. 

As expected, the results of vaccinating in the inter- 
val from 7 to 8 yr of age are disappointing if compared 
to the results of vaccinating from 1 to 2 yr of age, and 
vaccinating between 14 and 15 yr is almost useless. 
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Vaccination programmes against rubella were imple- 
mented in many countries (e.g., Refs. @, [HI El, HI, 
H3, However, vaccination coverages and strategies 

sometimes changed from one period to another. As men- 
tioned by Ukkonen [13], UK (1970) and Finland (1975) 
chose selective vaccination of 11- and 13-year-old girls to 
prevent rubella and such a strategy was not effective in 
eradicating the virus. These observed results agree with 
our simulation for vaccination from 14 to 15 yr of age. 
In 1998, rubella vaccine was introduced in Mexico into 
the childhood vaccination schedule at age 1 and 6 yr [2l[ , 
resulting in an intense decrease in the rubella incidence, 
in agreement with our simulations for vaccination from 1 
to 2 yr of age. 

B. Temporal evolution 

The simulations for the temporal evolution of the force 
of infection were based on the numerical solutions of 
the integral equation for A(a,i), using the parameters 
of P(a,a') for the Caieiras community 

Our first simulation considered a completely suscepti- 
ble population (which is not the case with Caieiras) . We 
then assumed that, at time t = 0, a proportion pi = 10~ 5 
of individuals with ages between 40 and 45 yr suddenly 
become infected. The resulting dynamics of the disease 
is shown in Fig. [SJ 

It can be noted that after a few oscillations the force 
of infection tends to the function Ao(a). We can also see 
that from around 40 yr onwards the force of infection 
stabilizes. Figure [7] displays a profile cut at 8 yr old. 

For our next simulation, the same conditions as the 
above simulation were applied, and a vaccination routine 
of form (|4ip with qq = 1 yr, a\ = 2 yr, to — 40 yr, and a 
vaccination coverage of 70% was added. The results for 
all ages are shown in Fig. [H It can be noted that after 
the introduction of the vaccination, the force of infection 
oscillates before reaching a steady state, much lower than 
Ao(o). The whole process takes around 40 yr to reach 
the new steady state. The fact that the process takes 
so long to reach a steady state does not recommend this 
vaccination strategy for eradicating the disease. 

The next simulation uses the same vaccination scheme, 
but with a vaccination coverage of 80%. The results for 
all ages are shown in Fig. [9] 

For this coverage one can see that the disease is erad- 
icated in approximately 20 yr. Again the fact that the 
process takes so long does not recommend this vaccina- 
tion strategy for control. 

C. Comparison of specific features with real data 

The strategy given by Eq. (|41~j) was actually adopted in 
the United States in 1969 [l(. The results of the impact 
on the incidence is shown in Fig. [10] (left-hand scale) 
together with the results of our simulation (right-hand 




FIG. 6: (a) Simulations for X(a,t), without vaccination, con- 
sidering a completely susceptible population, at time t — 0, 
except for a proportion 10~ 5 of individuals with ages between 
40 and 45 yr that suddenly become infected; (b) the same as 
(a), but with the time scale starting at t = 5 yr, so that the 
initial peak is cut and the resulting steady state is observable. 

scale) for a v that results in a 80% coverage. The model 
estimates for the number of new infections per 100 000 
population were calculated according to the following 
equation: 

y(t) = l^ da X(a,t)S(a,t) , (52) 

where J °° da X(a, t)S(a, t) is the number of new cases per 
unit time at time t. In the calculations, the upper limit 
of the above integral was taken to be 60 yr. 

It should be noted that the incidence calculated from 
the seroprevalence data is two orders of magnitude larger 
than the incidence that results from notification. In fact, 
it is known that only a fraction of all infections display 
the clinical features of rubella disease. In addition, only a 
fraction of those rubella cases is officially notified. How- 
ever, several qualitative features of the data are quite 
similar to those observed in the simulations described in 
the preceding section. Let us comment, in more detail, 
on the more significant similarities. 

In 1977, that is, 8 yr after the introduction of the pro- 
gram, it was noted that although the program was hav- 
ing a major impact on rubella in children, rubella rates 
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FIG. 7: Profile of Fig. [H cut at 8 yr age. The initial outbreak 
peak between t = and t = 3 yr almost exhausts the suscepti- 
ble fraction of the population. It takes about 3.5 years for the 
number of new susceptibles to accumulate in sufficient num- 
ber to trigger a second outbreak which eventually stabilizes 
at an endemic steady state. 



in those older than 15 yr were not substantially different 
from prevaccination rates. We shall see now that this 
effect is shown in our simulations. 

Figure ITT1 represents four cuts of Fig. [5J corresponding 
to 70% coverage, at the ages of 8, 16, 25 and 35 yr. It 
can be seen that the above mentioned effect is clearly 
observed. The drop in the force of infection at the age of 
8 yr is much larger than at 16 and 25 yr, and the effect 
at the age of 35 yr is almost negligible. 

However, about 15 years after the introduction of the 
vaccine, three major outbreaks are observed in the sim- 
ulations, and this may be dangerous. The pattern of 
several oscillations in the incidence of an infectious dis- 
ease, after the introduction of vaccination, has already 
been observed in real data [23| . If the vaccinal coverage 
is increased to 80%, after small outbreaks, the disease 
disappears, as shown in Fig. [9] 



V. SUMMARY 

In this paper, we analyzed the temporal evolution of 
the age-dependent force of infection and incidence of 
rubella, after the introduction of a very specific vaccina- 
tion program in a previously nonvaccinated population 
where rubella was in an endemic steady state. This very 
specific vaccination program consists in vaccinating chil- 
dren within a certain age range with a rate determined 
essentially by the public response to government adver- 



FIG. 8: Simulations for X(a,t), with the same initial condi- 
tions considered in Fig. [6] but including a vaccination routine 
of form (|41|) . with ao = 1 yr, ai = 2 yr, to = 40 yr, and a 
vaccination coverage of 70%. 



tisements. 

We conclude that a simple vaccination program is not 
very efficient (very slow) in the goal of eradicating the 
disease. This gives support to a mixed strategy proposed 
by Massad et al. Q, accepted and implemented by the 
government of the State of Sao Paulo. This strategy rec- 




FIG. 9: Vaccination scheme considered in the simulations of 
Fig. [S] but with a vaccination coverage of 80%. 
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FIG. 10: Impact on the number of reported cases (left-hand 
scale) of the vaccination strategy against rubella adopted in 
the U.S. in 1969 (data taken from CDC [2^1. together with 
the results of our simulation (right-hand scale) for a vaccina- 
tion rate that results in a 80% vaccination coverage. 
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ommended a mass vaccination campaign against rubella 
in the State of Sao Paulo for all children with ages be- 
tween 1 and 10 yr as an initial intervention followed by 
a vaccination program of the form given by (|41[) . in the 
routine calendar at 15 months of age. As reported in 
Refs. [sL fl3j| . the results were very good, and there was a 
considerable reduction in the number of rubella and con- 
genital rubella syndrome cases. The incidence of rubella 
and CRS remained at low levels with the routine vacci- 
nation program, in agreement with our simulation results 
for high vaccination coverages. 

We have also applied a formalism developed elsewhere 
Q to calculate the effects of vaccination routines de- 
signed to reduce or eliminate rubella. 

This formalism provides an integral equation for the 
force of infection in a steady state given the pattern of 
contacts between the members of the population and the 
specific form of the vaccination routines. 

To apply the formalism, the pattern of contacts be- 
tween the members of the population, the so-called con- 
tact function /3(a, a'), has to be estimated. Some sym- 
metries obeyed by (3{a, a') and a general form for it were 
studied in Sec. Ill Bl In Sec. 1111} the force of infection in 
the absence of the vaccination was calculated from sero- 
prevalence data from four communities in Brazil, Mex- 
ico, Finland and the United Kingdom. With this force 
of infection, in the absence of vaccination, the contact 
function f3(a, a') for each community was estimated. It 
was noted that f3(a, a') differed considerably between the 
communities studied, which is in agreement with the dif- 
ferences in the force of infection in the absence of vacci- 
nation. 

Finally, in Sec. IIV A[ the effects of several vaccina- 
tion routines were calculated for the four communities 
studied. As a general conclusion, one can say that vac- 
cination between 1 and 2 yr presents distinct advantages 
over any other strategy considered. In Caieiras, vacci- 
nation between 7 and 8 yr has the apparent advantage 
of shifting the average age of the first infection leftwards. 
However, if the coverage is above 60%, the impact of vac- 
cinating between 1 and 2 yr on the force of infection is 
twice as high as vaccinating between 7 and 8 yr. This re- 
sult confirms our previous analysis and recommendations 
of 1992 (Massad et al. [2|). In all other communities 
studied, vaccination between 7 and 8 yr results in very 
disappointing impact when compared with vaccination 
between 1 and 2 yr. 
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